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ABSTRACT 

Statistical properties of the intensity in adaptive optics images are usually modeled with a 
Rician distribution. We study the central point of the image, where this model is inappropriate 
for high to very high correction levels. The central point is an important problem because it 
gives the Strehl ratio distribution. We show that the central point distribution can be modeled 
using a non-central Y distribution. 

Subject headings: instrumentation: adaptive optics, instrumentation: high angular resolution 
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1. Introduction 

In this Letter, we study the statistics of the light intensity at the central point of adaptive optics (AO) 
images in the case of a high to very high A O correction. This pro blem corresponds to the statistical properties 
of the "instantaneous" Strehl Ratio (SR) (|Fusco fc Conanl 2004 ). The regime that we consider is relevant to 
current or future AO systems. 

Outside the central point, the statistical properties of the light intensity constrain the detection limits of 
faint companio ns to nearby stars. It can be described by a Rician distribution, bo th in the case of a partia l 
AO correction ( Canales fc Cagigal 1999a ). and in the case of a high or very high (lAime fc Soummer 2004 ). 



This model has been successfully v erified in Lick A O data hy brid j OrahaJ fcmd. and the c~of 
coronagraphic images is detailed in I Soummer et al.1 (|2007l ). 



T he intensity distribution at the center of AO-corrected images was previously studied bv lCagigal fc Canales 
(|l998l ): ICanales fc Cagigal (|200lh at very low correction levels, for SR i n the 10% — 20% regime. They showed 
a good agree ment between the Rician m odel, numerical simulations ( Canales fc Cagigal 1999bl ). and exper- 
imental data (jCagigal fc Canales! 120011 ). At high to very high AO correction levels, howev er, the Rician 



distri bution fails to model successfully the statistics of the intensity at the central point. Gladvsz et al 



(2006) provided evidence of this pro blem in short exposur e Lick AO images, and we further verify this issue 
using numerical simulations (Sec|3]). Idadvsz et al.l (|2006l ) modeled the observed SR distribution using the 
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statistical properties of the atmospheric seeing. In this paper, we describe an alternative approach based 
solely on image formation and considering stationary phase properties. The model con sists of a non-ce ntral 
r distribution for th e opposite of the modulus of the complex amplitude. Recently, iGladvsd ((2006) and 



Christ ou et al.l (|2006l ) independently developed similar approaches to modeling stationary processes, based 



on a model of the variance of the phase. 



2. Statistical properties of the intensity 



We briefly recall the derivatio n of the statistical properti es of AO images. A detailed presentation is 
given by lAime fc Soummerl (120041) an d ISoummer et al.l (|2007l ) , based on results known in the context of 



holography by Goodmanl (| 19751 . 2006 ). In this section we use a one-dimensional formalism for clarity, but 



the results are valid in the general case. The complex amplitude in the focal plane is the Fourier transform 
of the pupil plane complex amplitude: 



*(r) = J P(x)(A + a(x))t 



-2%7tx r 



dx, 



(i) 



where P(x) denotes the pupil function, normalized such as f P(x)dx = 1. A corresponds to the deterministic 
perfect part of the wave, and a(x) is the zero- mean complex error term: A J r a(x) = e l ^ x \ with A = E[e^^], 
where E denotes the expectation value. <p(x) denotes the phase of the wavefront in the pupil plane, assuming 
a zero- mean Gaussian, with a constant variance cr^. We use the classica l definition of SR as the ratio 
between the central point intensities in the actual and ideal case ( Hardy 1998). Using the extended Marechal 
approximation (SR « e _cr< ^), A becomes A = e -cr< ^ 2 , which is approximately the square root of SR. 

Outside the central point of the image, the distribution of the complex amplitude can be approximated 
using known results from signal processing. Assuming that the complex amplitude a(x) can be represented 
by discrete values, and that the correlation between two points in the pupil plane decreases with their relative 
distance, the probability density function (PDF) of th e complex ampli tude in the focal plane for the random 
part of EqQ] is an asymptotically circular Gaussian (|Brillingerlll98ll ). We recall that for a scalar circular 
Gaussian distribution, i.e. z ~ J\f c (0,o- 2 ), the real and imaginary parts are independent and have the same 
variance. We will consider that for r ^ 0, f follows a decentered Gaussian distribution \I/ ~ AT C (C, J s ), where 
C is the complex ampl itude in the focal plane corresponding to the perfect part of the wave A (see Fig. 2 of 
Aime fc Soummer (|200J)). The corresponding intensity follows a Rician distribution: 



Pi(i) 



1 



■ exp 



i + I c 



h 



'2Vn/7c 



i > 0, 



(2) 



where I c = |C| 2 , and Iq denotes the zero-order modified Bessel function of the first k ind. This model has 
been verifie d outside the central point of the image in simulations ( Soummer et al.ll2007 ) and in real adaptive 
optics data ([Fitzgerald fc GrahamI 20061 ) . 



The central point of the image is a particular case, which gives a natural estimato r of the SR. We tested 
the Rician distribution at this location with a numerical simulation, using PAOLA (jJolissaint et al.ll2006l ) 
to generate independent instantaneous intensity values. We chose the parameters of the AO systems to 
generate several sets of data with SR ranging between 50% and 95%. We used a maximum likelihood (ML) 
estimation of the parameters I c and 7 S , assuming Eq|2l The Likelihood is computed for the unbinned data 
and maximized using optimization routines of Mathematica. The starting parameters are obtained from 
the moments method. We performed both the y 2 test and the Kolomogorov-Smirnov (KS) test on these 
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results, using 10 identical bins for the x 2 tes t, and a Mont e- Carlo estimation of the KS distribution, since 
the parameters are estimated from the data ([jenkins 20031 ). Both tests conclude that the Rician model is 
incompatible with the data at the central point, with respective right tail values of 3.8 10 -4 and 1.6 10 -2 . 
This is confirmed by reproducing these tests for a few independent sets of simulated data. We illustrate the 
histogram of the simulation and the best fit obtained with the Rician model in FigQ3 This result is not 




Fig. 1. — Histogram of the central point intensity ( SR) for a numerical simulation of 1000 independent 
PAOLA phase screens with a SR of 85%. The error bars assume statistical Poisson noise. The best fit of 
a Rician PDF is obtained for I c = 0.83, I s = 3 10 -5 and is superimposed to the Histrogram. The Rician 
model fails to describe these data according to the x 2 and Kolomogorov- Smirnoff tests. 

surprising since the skewness of the Rician distribution is (2 + 3A)(1 + A) _3//2 , which is positive, whereas 
the data histogram cl early exhibits a negative skewness. Observations based on real AO data, carried out 
at Lick Observatory ( Gladvsz et al. 20061 ) and Palomar Observatory (R. Soummer & J. P. Lloyd 2007, in 



preparation), also reported negatively skewed SR distributions. 



3. Central point: Strehl Ratio distribution 

At the center of the image, the Fourier phase term of EqO] vanishes and the complex amplitude is simply 
the integral of the pupil complex amplitude, \I/(0) = J P(x) e^^dx. At low correction levels, the phase 4>{x) 
is large, and the vectors e l ^ x ^ take any orientation in the complex plane. The sum of a large number of these 
vectors over the aperture produces a random walk, and the corresponding complex amplitude distribution is 
asymptotically Gaussian, according to a central limit theorem. At high correction levels, the phase <j){x) is 
small and the vectors e l ^ x ^ are not oriented randomly in the complex plane. Their sum does not produce a 
random walk, and the corresponding distribution is not circular Gaussian. This is the case that we study in 
this section. When moving outside of the center, however, the Fourier phasors rotate the vectors e l ^ x ^ in the 
complex plane, and the sum of the resulting vectors produces a random walk. This explains qualitatively why 
we obtain a circular Gaussian distribution outside the central point, for r > A/D, even at very high correction 
levels. Figj2] shows the histogram of independent realizations of the complex amplitude, at the central point 
where it is clearly not circular Gaussian, and in the transition region (r <C \/D) where the circularization 
occurs. This phen omenon can also b e understood qualitatively by considering a Taylor expansion of the 
phase term e 1 ^ ([Perrin et all l2003h . At a very high SR, the expansion is mainly dominated by the first 



oder term i<j)(x) and the resulting distribution is an approximately imaginary Gaussian (see Figj2j left most 
panel). Fig |3] illustrates the complex amplitude histogram for a SR of about 80% where the effect of higher 
order terms is clearly visible, when compared to the left panel of Fig|2l At a high or very high SR, we can 
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Fig. 2. — Histogram of independent realizations of the complex amplitude in the focal plane for a SR of 
approximately 95% at three locations indicated in each panel. Left: The complex amplitude at the central 
point is not a circular Gaussian distribution. When moving away from the center, the distribution becomes 
progressively circular, and is fully circularized for positions r > X/D (not represented here). 
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Fig. 3. — Histogram of independent realizations of the complex amplitude at the central point in the focal 
plane, for a SR of 80%. The corresponding distribution of intensity is negatively skewed. This can be seen 
by looking closely at the distribution of the modulus, along the real axis for example. 
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limit the expansion of e l ^ x ^ to the second term ( Sivaramakrishnan et al. 2QQ2h and we obtain the intensity: 



1 



/ P{x) 



^xf 



dx 



P{x)(j){x)dx 



(3) 



The second term corresponds to the mean of the phase residuals over the pupil. It seems reasonable to neglect 
this term compared to the first term, assuming the phase to be zero mea JE Given that 1 > J P{x)^-dx, 
we can revert to the complex amplitude: 



*(0) 



1 



P(x) 



<j>{xf 



dx 



(4) 



by introducing a random phase term cp. Note that ip is different from the pupil phase 



The square of the phase (f)(x) 2 is proportional to a xi distributed random variable: (j)(x) 2 = a 2 Q{x) 
where Q(x) ~ Xi- The integral U = J P(x)(f)(x) 2 dx can be reasonably approximated by the discrete sum: 
U = <j' 2 Yli=i Q( x i)i where the number of discrete samples can be assumed as large as necessary, and where 
a' 2 = Sa 2 with S the uniform integration step. The variable Yli=i Q( x i) ~ x|> where k denotes the degree 
of freedom of the x 2 distribution (he re corresponding to the number of discrete samples). The distribution 
of U is a T-distribution T(k/2, 2a' 2 ) \ Johnson et al.lll994h : 



Pu(u) 



(2a f2 )- k / 2 u k / 2 - 1 e- u ^ 2 
W2) 



u > 0. 



(5) 



The opposite of the modulus of the complex amplitude A = U/2 — 1 defined as ^(0) = — Ae Jlf thus follows 
a three-parameter non-central Gamma distribution r(fc/2, <j /2 , — 1). The approximation of EqJH sets an 
arbitrary constraint on the m ean of the dis tribution: E(|\I/(0)|) = 1 — cr^/2, which corresponds to the regime 
of Marechal's approximation (|Hardyl ll998). The constraint on the mean can be released without changing 
the statistical model by relaxing the third parameter. Finally, we model A by the T(a,/3, 7) distribution 
defined as: 

(a-7)-- 1 exp(-(a- 7 )//3) 



The mean and variance of this distribution are: 



-, a > 0, (3 > 0, a > 7. 



(6) 



E[A] =a(3 + <y, var[^4] = a(3 2 . 

The skewness of the distribution is 2/y/a which is i ndeed negative for \f l 
parameters can be obtained by a moments method ( Johnson et al. 1994h : 



(7) 



-A. Estimations of the three 



4^=^,7 = Mi- M 



Ms 



2fi 



Ms 



(8) 



where /ifc is the classical empirical estimator of the kth centered moment. I Johnson et al.l (|1994l ) gives an 
algorithm for the computation of the ML estimation of a, /?, and 7 together with their asymptotic variance. 
The SR can be expressed directly using the mean and variance of the distribution (EqEj): 



SR=E[A 2 ] = a(3 2 + M + 7) 2 



(9) 



x Note that whereas this term can be neglected in the intensity, it cannot be neglected in the second order approximation of 
the complex amplitude, because it corresponds to the imaginary part, which cannot be neglected compared to the real part. 
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We obtain a new estimator of the SR by substituting the parameters in Eq|9]with their estimated values, using 
the moments method or the ML. Using the moments method, the estimation of the SR is: SR = /if + 2/i|. 
FigHl shows the results of numerical simulations for SR « 90%, including 2000 independent realizations of 
phase and amplitude screens, corresponding to the atmospheric corrected wavefronts and to the atmospheric 
scintillation. We perform a ML estimation of the three parameters a, /3, and 7, using the moments method 
to set the initial values. The estimated SR (89, 58%) from Eqj9] matches the simulation parameters very well 
(S# = 89.50%). 



200 




Fig. 4. — Numerical simulation for SR~ 90%. The simulation include both amplitude and phase screens. 
The figure shows the histogram of \fl and the proposed model. Both methods (moments, and maximum 
likelihood) give almost identical results. 

We can derive an analytical expression for the PDF of the intensity I = A 2 , using pa{°) from Eq|6j 

Pi(i) = ^= (pa(VI) +PA{-Vi)) ,i > (10) 

Fig|5] shows two typical PDF for the intensity, where the parameters have been estimated from simulated 
data, using the Gamma distribution for the modulus. 
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Fig. 5. — Examples of two analytical PDF models for the intensity for estimated SR of 89, 58% and 94.50%. 
The skewness of the distribution decreases with SR and is almost zero at very high SR. 
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4. Conclusion 

The Rician distribution, which describes successfully the statistics of AO images outside of the central 
point, is not appropriate for the central point of the image at high correction levels. This specific location is 
particularly important as it corresponds to the SR. At high SR levels, the complex amplitude at the central 
point is not described by a circular Gaussian distribution, which can be explained by the role of Fourier 
phasors that vanish at the central point and circularize the distribution outside the center. Numerical 
simulations support this conclusion. 

We propose an alternative approach, based on a model of the modulus of the complex amplitude, well 
described by a non-central Gamma distribution. Our method enables the use of standard algorithms for the 
estimation of the parameters of this well known distribution. When the parameters have been estimated 
from the PDF of the modulus, it is possible to revert to the PDF for the intensity by the appropriate 
transformation. Application to real data sets is und er study. It would be interesting to compar e the approach 
based on the seeing statistics (jdadvsz et al.ll2006l ) or on the model of the phase variance ( Christou et al 



20061 ). with our results which do not require fluctuations of the seeing. 
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